Downloaded from https://www.pnas.org by 178.204.250.13 on January 11, 2023 from IP address 178.204.250.13. 


PNAS 


RESEARCH ARTICLE 


EARTH, ATMOSPHERIC, AND PLANETARY SCIENCES 


af’ OPEN ACCESS ® 


Check for | 
| updates 


The Bering Strait was flooded 10,000 years before the Last 


Glacial Maximum 


Jesse R. Farmer?" © Tamara Pico®' © Ona M. Underwood’, Rebecca Cleveland Stout??®, Julie Granger ®, Thomas M. Cronin’®, 


Francois Fripiat®?®, Alfredo Martinez-Garcia® ®, Gerald H. Haug?" ®, and Daniel M. Sigman? ® 


Edited by Claire Waelbroeck, Laboratoire d'Océanographie et du Climat: Expérimentations et Approches Numériques, Institut Pierre-Simon Laplace, 


Paris, France; received April 18, 2022; accepted October 21, 2022 by Editorial Board Member Jean Jouzel 


‘The cyclic growth and decay of continental ice sheets can be reconstructed from the 
history of global sea level. Sea level is relatively well constrained for the Last Glacial 
Maximum (LGM, 26,500 to 19,000 y ago, 26.5 to 19 ka) and the ensuing deglaciation. 
However, sea-level estimates for the period of ice-sheet growth before the LGM vary 
by > 60 m, an uncertainty comparable to the sea-level equivalent of the contemporary 
Antarctic Ice Sheet. Here, we constrain sea level prior to the LGM by reconstructing the 
flooding history of the shallow Bering Strait since 46 ka. Using a geochemical proxy of 
Pacific nutrient input to the Arctic Ocean, we find that the Bering Strait was flooded 
from the beginning of our records at 46 ka until 35.714 ka. To match this flooding 


history, our sea-level model requires an ice history in which over 50% of the LGM’s 
global peak ice volume grew after 46 ka. This finding implies that global ice volume 


and climate were not linearly coupled during the last ice age, with implications for the 
controls on each. Moreover, our results shorten the time window between the opening 
of the Bering Land Bridge and the arrival of humans in the Americas. 


Arctic Ocean | Bering Strait | sea level | foraminifera-bound N isotopes | glacial isostatic 
adjustment 


The Bering Strait, the -53 m deep ocean passage that separates Asia from North America 
(1), is the only Northern Hemisphere connection between the Pacific and Atlantic Oceans. 
Today, -1 Sv (10° m° s`!) of low-salinity Pacific seawater flows northward across the Bering 
Strait (2) and contributes to the relative freshness of the upper Arctic Ocean (Fig. 1). 
Export of these Arctic waters into the subpolar North Atlantic modifies the surface waters 
that form North Atlantic Deep Water (3) and may cause feedbacks between the North 
Pacific and North Atlantic Oceans (4—6). 

During the last glacial cycle, sea-level changes driven by the growth and decay of con- 
tinental ice sheets exposed and flooded the Bering Strait. The Bering Strait flooded most 
recently between 13 and 11 ka during sea-level rise caused by the melting of ice sheets 
(11, 12). Before 13 ka, the Bering Strait was subaerially exposed due to lowered sea level 
from extensive continental ice sheets during the Last Glacial Maximum (LGM) (26.5 to 
19 ka, ref. (13)) and early deglaciation. At this time, Asia and North America were con- 
nected by the Bering Land Bridge, a proposed route by which human populations first 
entered the Americas (14-16). However, there is great uncertainty as to the timing with 
which ice growth leading up to the LGM exposed the Bering Strait and formed the Bering 
Land Bridge. Hopkins (14) initially reported geological evidence for Bering Strait sub- 
mergence prior to the LGM, but this evidence is debated (17). 

Relative sea level at the Bering Strait (RSL,,<,) is affected by global mean sea level (GMSL) 
and the solid Earth response to the growth and decay of ice sheets through the process of 
glacial isostatic adjustment (11). Uncertainty regarding Bering Strait submergence prior to 
the LGM reflects the correspondingly high uncertainty as to the history of GMSL leading 
up to the LGM. Estimates of GMSL between 50 and 30 ka from various geological and 
geochemical data range between -25 and -105 m (Fig. 2B, refs. (18-24)). This GMSL 
uncertainty of > 60 m, which exceeds the entire sea-level equivalent of the modern Antarctic 
Ice Sheet (58 m, (25)), reflects in part the paucity of geological sea-level observations: 
Advancing ice sheets razed evidence of prior ice margins (26), and sea-level rise during the 
last deglaciation destroyed or submerged ancient coastlines. In addition, radiometrically 
dated coral sea-level markers are limited to uplifted terraces during this time, with recon- 
structed paleo-elevations that are subject to errors in uplift corrections (27, 28). 

The submergence history of the Bering Strait is important in diverse contexts. First, 
considering its modern sill depth of -53 m, the Bering Strait’s submergence history, when 
corrected for glacial isostatic adjustment, could serve as a critically needed constraint on 
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GMSL between 50 and 30 ka. IF RSLges is closely correlated with 
GMSL, sea-level inferences from the ICE-5G global ice sheet 
history (19) (black line, Fig. 2B), oxygen isotope records in the 
Red Sea (20) (red line, Fig. 2B), and a stacked sea-level equivalent 
oxygen isotope record (21) (blue line, Fig. 2B) suggest that the 
Bering Strait would have been subaerially exposed during this 
time, whereas recent GMSL estimates based on glacial isostatic 
adjustment analyses of sea-level data (22, 23) (purple line, Fig. 2B) 
and ice margin constraints (24) (gray line, Fig. 2B) indicate that 
the Bering Strait would have been submerged. Second, the Bering 
Strait’s sea-level history impacts the connectivity of the Pacific and 
Atlantic Oceans, which has been hypothesized to control abrupt 
(millennial-scale) climate and ocean circulation variability during 
the last ice age (5, 6). Third, the submergence history affects when 
terrestrial migrations were possible between Asia and North 
America across the Bering Land Bridge, putatively leading to the 
first arrival of humans in North America. 

Here, we reconstruct the history of Bering Strait submergence 
since 46 ka using a geochemical proxy for Pacific water input to 
the western Arctic Ocean (8) and simulations of relative sea level 
at the Bering Strait (11, 30). These reconstructions show that, 
contrary to previous assumptions (5, 6) and congruent with the 
original hypothesis of Hopkins (14), the Bering Strait was flooded 
by at least 46 ka and that the Bering Land Bridge formed only 
after -36 ka. 

The geochemical proxy for Bering Strait submergence is based 
on regional features of the marine nitrogen (N) cycle in the polar 
Northern hemisphere oceans. The N isotopic composition (6 ‘N= 
ICN N) ample N/N) — 1] * 1000) of nitrate (NO3), the 
primary form of fixed N supplied to surface ocean ecosystems, 
varies by about 3 parts per thousand (3%) among the high-lati- 
tude North Atlantic, North Pacific, and western Arctic Oceans 
today (Fig. 1A). Specifically, the 6'°N of nitrate supplied to the 
surface mixed layer is substantially higher in the western Arctic 
(~8%o) than in the eastern Arctic or subpolar North Atlantic 
(~5%o). The high western Arctic nitrate N arises from two 
features unique to the Pacific inflow across the Bering Strait. First, 
the nitrate inflow across the Bering Strait is elevated in 5N rel- 
ative to the high-latitude North Atlantic Ocean. This reflects both 
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the higher nitrate 5°N in subarctic North Pacific subsurface 
waters relative to subarctic North Atlantic waters (SI Appendix, 
Fig. S1) due to water column denitrification in the North Pacific 
interior, and partial nitrate assimilation in Bering Sea surface 
waters (32). Second, the high nitrate concentration of this Bering 
Strait inflow fuels high shelf productivity, which in turn initiates 
coupled partial nitrification-denitrification (CPND) on the Bering 
Sea shelf and the western Arctic shelves; this CPND further ele- 
vates nitrate 8 N (33-36). Importantly, the Arctic Ocean’s CPND 
is sited exclusively in areas that are influenced by the Bering Strait 
inflow (33-36) (SI Appendix, section S2 and Figs. S1 and S2). 
Accordingly, the subsurface nitrate 6N elevation in the western 
Arctic is directly tied to Pacific nitrate input (Fig. LA). 

We reconstruct the N isotope signature associated with the 
Bering Strait inflow back to 46 ka by measuring the ôN of organic 
matter bound within the planktonic foraminifer Neogloboquadrina 
pachyderma (hereafter, 89N y p.) from three sediment cores in the 
western Arctic Ocean and, as a control, from one core in the central 
Arctic Ocean outside the direct influence of the Bering Strait inflow 
(Fig. LA). Planktonic foraminifera-bound 5°N (such as SÝN Np) 
reflects the 6’°N of organic matter produced in surface waters 
(8, 37), which depends on the 6’°N of the subsurface nitrate supply 
and the degree of nitrate consumption, that is, the summertime 
drawdown of nitrate as a proportion of annual nitrate supply (38). 
Summertime nitrate consumption is complete in the western Arctic 
today due to the highly stratified upper water column (39), which 
arises in part from continual supply of low-salinity seawater across 
the Bering Strait (8) (Fig. 1). 

Farmer et al. (8) showed that western Arctic 6'°N y, p. recorded 
the most recent postglacial flooding of the Bering Strait, which 
has been dated by independent methods to 13 to 11 ka (e.g., refs. 
(11) and (12)) (Fig. 2D). During Marine Isotope Stage (MIS) 2 
(29 to 11.7 ka, including the LGM), western Arctic SPN yy 
ranged from 4.5 to 6%o, requiring Atlantic-sourced nitrate as well 
as incomplete nitrate consumption due to weaker density strati- 
fication (8). Around 11.5 ka, western Arctic Oo Nis rapidly 
rose to values of 7.8 to 8.8%o throughout the Holocene. This 
OO Nis rise resulted from flooding of the Bering Strait, which 
introduced nitrate-rich, high-6 'N Pacific waters, triggered shelf 


Fig. 1. Hydrography of the polar Northern Hemisphere and Bering Strait bathymetry. (A) Mixed layer (10 m depth) salinity (7) (blue shading) and subsurface 
(~50 to 200 m) nitrogen isotopic composition (5'°N) of nitrate (diamonds). From published data sets, sample depths were chosen to capture the nitrate being 
supplied to the mixed layer in the Arctic and North Atlantic (NA), and the North Pacific (NP) nitrate being transported across the Bering Strait (see data sources 
and selection criteria in S/ Appendix, section S2). The mean salinity and nitrate 5'°N values for the western Arctic (WA), NP, and NA are indicated on the color bars. 
Cyan circles show locations of sediment cores; light blue and orange arrows show schematic circulation of Atlantic- and Pacific-sourced nitrate-rich subsurface 
waters, respectively (8, 9). (B) Bering Strait bathymetry (1) contoured at 10 m intervals. Red lines indicate principal ocean transport pathways (10). 


2of7 https://doi.org/10.1073/pnas.2206742119 


pnas.org 


Downloaded from https://www.pnas.org by 178.204.250.13 on January 11, 2023 from IP address 178.204.250.13. 


A ° 10 20 30 40 50 
-30 
= MIS 1 MIS 2 MIS 3 
X 
O -35 
S 
© 
a 
© -40 
[0] 
m 
-45 
B 2 
e 
5 Fossil corals A 
Moa a Huon (corr.) A 2 
a 
— PaleoMIST O 
— Red Sea z 
— PCA stack D 
C — ICE-5G i? 
= 3 
S -120 
(si 
s 
ke] 
2 
E 
= 
E 
© 
a 


ČN ate “sa 0%) “Neg Ono 


Age (ka) 


Fig. 2. Arctic Ocean foraminifera-bound N isotope records compared with 
related records of climate and sea-level change since 50 ka. (A) NGRIP ice core 
5'°0, reflecting Greenland air temperature (29). (B) Sea-level reconstructions 
from ICE-5G (19) (black), the Red Sea (20) (blue-red shading), stacked oxygen 
isotope-based sea-level reconstructions (21) (purple with purple shading), 
PaleoMIST (24) (orange, dashed line is minimal MIS 3 GMSL scenario), ICE- 
PC (red, this study and ref. (30)), previously published MIS 3 coral sea-level 
benchmarks (27) (gray triangles), and recently corrected Huon Peninsula sea- 
level datums (28) (yellow triangles). (C) Global benthic foraminifera 8'80 stack 
(31), which records both global ice volume and deep ocean temperature and 
is widely applied as a sea-level proxy. (D) Western (circles) and central (gray 
diamonds) Arctic Ocean 5'°Ny,, (this study); vertical error bars denote the 
larger of measured or long-term replicate 5’ Nyp, precision (Methods), and 
horizontal error bars denote 68% quantiles (equivalent to +1sd) of the age- 
depth model. Timing of MIS 1 to 3 according to (31) is denoted at the top. 


CPND that further elevated nitrate 6'°N in the western Arctic, 
and stratified the upper water column, leading to complete 
nitrate consumption (8). 

New 6? Ny p. records from three western Arctic Ocean sediment 
cores extend these reconstructions to the limit of radiocarbon 
dating (~50 ka), through the pre-LGM interval of larger sea-level 
uncertainty (Fig. 2). Results exhibit three distinct intervals sepa- 
rated by two rapid transitions. The youngest transition, a6'’’N jy, 
rise of ~3%o around 11 ka, records the postglacial flooding of the 
Bering Strait discussed above, ending the low 5'°N p. (of 4.5 to 
6%o) of MIS 2. The new data show that the low ô' Nw, of MIS 
2 did not extend back through MIS 3 (57 to 29 ka). Instead, all 
the three western Arctic sites show high 6'°N np. Values of 8 to 
10.5%o prior to 35 to 40 ka; these values are equivalent to or 
higher than Holocene 6N y p. at these locations. Moreover, west- 
ern Arctic 6 -Ny pvalues are 2 to 3.5%o higher than those meas- 
ured in the central Arctic before 35 to 40 ka (Fig. 2D). Western 
and central Arctic 6 -Ny, „Values converge after 35 ka, with similar 
6° Ny » Values (of 4.5 to 6%o) lasting in both regions until 11 ka. 

Moving forward in time from MIS 3, all the three western 
Arctic cores show a ~3%0 5 PN y p. decline within 2 cm of sediment 
(equating to 1 to 2 kyr based on age models; SI Appendix, section 
S1 and Fig. S3). Radiocarbon age models (SI Appendix, Fig. S3) 
date this transition to 40 ka at B8 and 35 ka at B12 and B17 
(SI Appendix, Fig. S5A). The data are consistent with a simultaneous 
6° Ny p decline in all the three cores, given the age model 
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uncertainties of up to +3 ka (95% CI) (SI Appendix, Fig. S5) and 
also considering the potential effects of bioturbation at these low 
sedimentation rates. Conversely, it would be difficult to explain 
different timings in the large and similar magnitude of 6 NK, 
decline at these three proximal western Arctic cores (which are all 
located within 360 km of one another; Fig. 14). Moreover, true 
diachrony in the 5°-Ny p. decline among sites is inconsistent with 
the decline in 6°-N y, p. occurring first at Site B8, which is the site 
closest to the Bering Strait (Fig. LA) and so would presumably be 
the last site to lose an isotopic signal emanating from the Strait. 
Thus, we consider that the 6N y, decline is contemporaneous 
at these locations and assign a median (+ interquartile range) age 
for the transition of 35.7133 ka (SI Appendix, Fig. S5A). 

The MIS 3 western Arctic 6 VN y p. records require the presence 
ofa high-6'°N nitrate source that was rapidly removed around 36 
ka. A terrigenous N source can be excluded, as both dissolved and 
particulate nitrogen inputs from Arctic rivers are too low in 5°N 
[2 to 5%o; (40, 41)] to explain the elevated western Arctic 5Nyy r 
values. This is also supported by low nitrate 8N (<5%o0) in the 
Kara Sea and the Laptev Sea, where large terrigenous N contribu- 
tions are expected (Fig. 1A and SI Appendix, Fig. S2). Additionally, 
our central Arctic 5'°N _ record (gray diamonds in Fig. 2D) 
requires that the high-6 IN nitrate was limited to the western 
Arctic during MIS 3, as occurs today (Fig. 1). The central Arctic 
6PNy p. record averages 7.1%o prior to 35 ka; this value is con- 
sistent with foraminifera-bound 6'°N records from the North 
Atlantic that indicate a regional upper water column nitrate 6 DN 
of 5 to 6%o during MIS 3 (42). Thus, the central Arctic most 
likely received nitrate from the North Atlantic Ocean during MIS 
3, as occurred throughout the last 35 kyr (8). The 6° Nx, „ homo- 
geneity of the western and central Arctic between 35 and 11 ka 
(Fig. 2D) points to similar North Atlantic nitrate sources and a 
shared condition of weak upper ocean stratification in the two 
regions (8). Finally, the spatial gradient in 6°N Np. between the 
western and central Arctic that occurred prior to ~36 ka redevel- 
oped by 11 ka with the deglacial flooding of the Bering Strait. 
Thus, from an N isotopic perspective, MIS 3 and the Holocene 
appear remarkably similar in the Arctic (Fig. 2). 

Given the above evidence, we conclude that the Bering Strait was 
flooded prior to -36 ka. With a flooded Bering Strait, high-6'°N 
nitrate from the Bering Sea would have been transported northward 
into the western Arctic Ocean as occurs today (33-36) (SI Appendix, 
Fig. S2). High nutrient concentrations in these Bering Strait inflow 
waters would also have fueled high primary productivity at the shelf 
break, triggering CPND that further elevates nitrate 6 DN in the 
western Arctic today (Fig. 1 and SI Appendix, section S3 and Fig. 
S1D). Finally, nitrate consumption in the western Arctic is complete 
due to the strong density stratification of the region, which appears 
contingent on the low-salinity Bering Strait inflow (8). The exist- 
ence of this inflow before 36 ka may have strengthened western 
Arctic stratification at that time, leading to more complete surface 
ocean nitrate consumption in the western Arctic, and this may be 
required to reach the high 6'°N np. of MIS 3. 

It is noteworthy that western Arctic 6'°N y p. is higher between 
40 and 46 ka than during the Holocene. This could reflect a greater 
extent of halocline and surface waters influenced by the Bering 
Strait inflow in the former period, enhanced CPND on the Bering 
Sea and western Arctic shelves, and/or a higher 5°N of the sub- 
arctic North Pacific nitrate source flowing onto the Bering Sea 
shelf due to more complete summertime nitrate consumption in 
the Bering Sea, as suggested by higher diatom-bound 6’’N in the 
Bering Sea during MIS 3 (43). While our current data do not 
distinguish among these explanations, a flooded Bering Strait 
during MIS 3 is required in all cases. 
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After ~36 ka, the rapid 8” Ny,, decline at all three western 
Arctic sites and 6 -N y p. Values of <7%o indicate the cessation of 
the Bering Strait inflow to the Arctic, following the same logic as 
outlined above (SI Appendix, section S3). This reflects the subaerial 
exposure of the Bering Strait and thus formation of the Bering 
Land Bridge, with a terrestrial connection persisting through the 
LGM (11, 12, 14, 17). 

To assess the quantitative implications of our findings for global 
ice volume during MIS 3, we model glacial isostatic adjustment 
and relative sea level at the Bering Strait (RSLpes). Our simulations 
assume that the modern sill depth of the Bering Strait (-53 m) has 
not substantially changed since 50 ka due to vertical displacement 
from processes unrelated to ice loading, such as longer term local 
tectonic, erosion, and sedimentation effects, as these are unlikely 
to be significant over the short time interval of our study (e.g., ref. 
(11)). We performed gravitationally self-consistent glacial isostatic 
adjustment simulations using three ice histories that encompass 
the range of MIS 3 GMSL estimates (Figs. 2B and 3B). The first 
ice history, ICE-5G (19), is characterized by MIS 3 GMSL values 
of -82 to -100 m (Fig. 3B). Peak predicted RSL,.<, during mid- 
MIS 3 (50 to 35 ka) is -65 m, below the modern sill depth (-53 
m, Fig. 3C). The second i ice history is derived from a GMSL history 
constructed by scaling 5'%O records to coral reef sea-level records 
((18); Materials and Methods) and is characterized by a GMSL of 
-71 to -62 m during mid-MIS 3. RSL,.<, based on this ice history 
(18) predicts an exposed Bering Strait from 50 to 35 ka, with 
RSLpes: below or equal to the sill depth (-58 to -53 m) (Fig. 3C). 
In summary, neither of these ice volume reconstructions associated 
with canonical sea-level histories (18, 19) are consistent with our 
finding of a flooded Bering Strait during MIS 3. 

Recent studies on the magnitude of GMSL during mid-MIS 
3 (50 to 35 ka) suggest substantially higher peak GMSL compared 
with the previous reconstructions (22—24). We estimate RSLgest 
for a set of ice histories consistent with these recent GMSL con- 
straints. One important caveat to these alternative ice histories is 
that RSLpes: is sensitive to the nearby ice sheet history (specifically, 
the Cordilleran Ice Sheet (CIS)) (11, 30). The CIS extent during 
MIS 3 is poorly constrained by field data, representing a source 
of uncertainty (44) (SI Appendix, section S4). To address this, we 
simulate three cases of MIS 3 CIS geometry (large, intermediate, 
small; SI Appendix, section S4), each of which maintains the same 
ICE-PC-derived GMSL history (30). In all three simulations, the 
Bering Strait is flooded during MIS 3 until after -34 ka 
(SI Appendix, Fig. S7). For the intermediate CIS history, RSL».<, 
is -11 to -20 m (yielding a water depth of 33 to 42 m; Fig. 3C) 
during mid-MIS 3. Our relative sea-level predictions suggest that 
a greater extent of the CIS results in a more deeply submerged 
Bering Strait (that is, a higher RSL».<,). Nevertheless, the three 
CIS geometries produce RSL».s, predictions that differ by less 
than 10 m (Materials and Methods and SI Appendix, Fig. S7). 
During MIS 3 and the LGM, glacial isostatic adjustment causes 
RSLpes: to be higher than GMSL by 20 to 30 m during MIS 3 
and the LGM (compare Fig. 3 B and C). Nevertheless, our sim- 
ulations show that GMSL changes are the dominant control on 
relative sea-level changes at the Bering Strait leading up to the 
LGM. Thus, a fully flooded Bering Strait during MIS 3 appears 
to require a global ice volume history in which MIS 3 ice volume 
is <50% (and in the case of ICE-PC, <30%) that of MIS 2 (22- 
24). This ice history contrasts with the long-dominant view of 
only modest ice retreat during MIS 3 (e.g., refs. (19-21)). A 
flooded Bering Strait before ~36 ka also requires that ice sheets 
grew quickly from MIS 3 until the LGM (Fig. 3B), consistent 
with previous findings (23, 45). 
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Fig. 3. Western Arctic Ocean 5° Nvps global mean sea-level reconstructions, 
and glacial isostatic adjustment simulation of Bering Strait relative sea-level 
history. (A) Western Arctic Ocean 5°Nyps ; dashed horizontal bar denotes 
maximum expected ô' Ny. for an exposed Bering Strait (Methods). Vertical 
error bars denote the larger of measured or long-term replicate 5' 5Nyp. 
precision; horizontal error bars denote 68% quantiles (equivalent to +1sd) of 
the age-depth model. (B) GMSL from ICE-5G (19) (black), from an ice history 
constructed from the pre-LGM GMSL in ref. (18) (W02-derived, blue), and from 
ICE-PC (red). (C) Relative sea level at the Bering Strait from ice volume histories 
in B, where ICE-PC RSLgeg is based on the intermediate Cordilleran Ice Sheet 
history (Materials and Methods). Black dashed line in B and C denotes modern 
sill depth of the Bering Strait (-53 m). Gray vertical shading denotes the timing 
of pre-LGM Bering Strait flooding reconstructed from panel (A). Colored bar 
at the bottom of C shows interpreted Bering Strait sea-level history. Brown 
vertical line denotes the reconstructed timing of Bering Strait closure in MIS 3; 
dashed brown lines are +95% confidence intervals (this study). Blue vertical line 
and dashes denote the mean timing and range of observations for postglacial 
Bering Strait flooding, respectively (8, 11, 12). 


Our Bering Strait submergence and GMSL reconstructions have 
implications for climate. First, relative sea level above -53 m at the 
Bering Strait until -36 ka (Fig. 44) supports recent estimates of peak 
MIS 3 GMSL near -40 m (22-24, 46, 47), indicating that peak 
MIS 3 global ice volumes were more similar to the Holocene than 
to the LGM (Fig. 48). In contrast, global temperature proxies sug- 
gest that MIS 3 was notably colder than the Holocene and only 
slightly warmer than the LGM (Fig. 4D). If global temperatures 
during MIS 3 were indeed similar to the LGM, why was ice volume 
so much lower during MIS 3? The growth of ice volume in late MIS 
3 might reflect nonlinearity in the sensitivity of ice volume to CO, 
radiative forcing (Fig. 4£), with ice volume responding more to CO, 
change when CO, is low and climate is cold (48). However, such 
nonlinearity would need to be very strong to explain the observa- 
tions, with the radiative forcing from the decline in atmospheric 
CO, from 220 to 190 ppmv between 40 and 30 ka driving >60 m 
sea-level equivalent ice volume growth (Fig. 4 Band £). An alterna- 
tive possibility is that Northern Hemisphere ice volume growth was 
more directly controlled by peak summer insolation than by global 
temperature (13, 24, 49), as supported by the coincidence of declin- 
ing 65°N summer insolation (Fig. 4C) and increasing ice volume 
through late MIS 3 (Fig. 4B) despite relatively little global cooling 
(Fig. 4D). This distinction may indicate a strong role for local sum- 
mer insolation in driving ice loss (Fig. 4 B and C), for example, by 
direct irradiance-driven ablation (e.g., ref. (50)) or as a local driver 
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Fig. 4. Comparison of histories of Bering Strait flooding, sea level, climate, 
and climate forcings over the last 50 kyr. (A) Bering Strait and Bering Land 
Bridge intervals based on western Arctic 5 Nwp, (as in Fig. 3); (B) ICE-PC (this 
study, red) and PaleoMIST ((24), orange) sea-level reconstructions. (C) 65°N 
summer solstice mean daily insolation (dashed, W m~?) and caloric summer 
half-year insolation (51) (solid, GJ m). (D) Global average surface temperature 
anomaly (52) (brown), Antarctic ice core temperature stack (ATS) anomaly 
(53) (cyan), and global sea surface temperature (SST) anomaly (49) (blue, 
shading is 95% Cl). (E) Atmospheric CO, concentration (54, 55) (dark red) and 
cumulative greenhouse gas (CO,, CH,, and N,O) radiative forcing (56) (ARFeyc 
black, shading is +1sd). 


of summertime air temperature change over the Northern 
Hemisphere ice sheets. 

Second, given the importance of ice albedo to global temperature, 
a complementary question emerges: If ice volume during MIS 3 
was sufficiently low to allow for a flooded Bering Strait, then why 
was MIS 3 so cold? One possibility is that the radiative impact of 
low ice volume was compensated by spatially extensive but thin 
land ice and snow cover as well as sea ice during MIS 3 (49). The 
cold MIS 3 climate could also indicate an outsized importance of 
low atmospheric CO, concentrations during MIS 3 (200 to 220 
ppmv, Fig. 4£) to Earth’s radiative balance (49). In general, the 
deviation from linear correlation between climate and ice volume 
during MIS 3 provides a new test for models of climate and glaci- 
ation. Extending the use of the Bering Strait as a sea-level gauge 
further back in time will likely offer additional constraints on the 
mechanisms underlying the glacial cycles as glaciological and cli- 
matic phenomena. 

Third, in both simple and fully coupled climate models, the 
Bering Strait’s sea-level history contributes to the response of North 
Atlantic Deep Water production to regional freshwater perturba- 
tions (5, 6). These simulations led to the proposal that subaerial 
exposure of the Bering Strait was a prerequisite for the millenni- 
al-scale Northern hemisphere climate variability during the last 
glacial period (5, 6). However, our data indicate that millennial-scale 
variability in climate and ocean circulation occurred both when the 


Bering Strait was flooded (46 to 36 ka) and exposed (36 to 11.5 ka) 
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(Fig. 2). Thus, closure of the Bering Strait was not needed for the 
occurrence of millennial-scale variability during the last ice age. 

Finally, our reconstructed Bering Strait submergence history 
has implications for humans' arrival in the Americas. Recent evi- 
dence for human presence in North America between 23 and 21 
ka (15) appears to require that humans migrated across Beringia 
before the onset of LGM conditions (15, 16). Limited regional 
archaeological records of human occupation before and during 
the LGM, particularly from eastern Beringia, preclude identifica- 
tion of the key factors controlling human migration at this time. 
However, comparison of our results with demographic age mod- 
eling of Beringian populations suggests one possibility. Humans 
were present in North Siberia as early as 45 ka (57), with a unique 
Ancient North Siberian population diverging from East Asians by 
39 ka (95% CI 45.8 to 32.2 ka, ref. (58)). A second distinct 
Ancient Beringian population with direct ancestry to Native 
Americans emerged from East Asians at 36 + 1.5 ka (16). Our 
data indicate that Siberia and Alaska were separated by the Bering 
Strait until the Bering Land Bridge formed at 7 ka. Together, 
these findings suggest that, in the context of the last -50 kyr, 
humans migrated into the Americas as soon as the Bering Land 
Bridge allowed for their passage. If correct, such simultaneity 
implies a strong drive for migration among the ice-age human 
populations in Siberia. 


Materials and Methods 


Three western Arctic cores were obtained from Mendeleev Ridge: Site B8 
(78.13°N, 176.74°W, 1,031 m water depth), Site B12 (79.99°N, 174.29°W, 
1,609 m water depth), and Site B17 (81.27°N, 178.97°E, 2,217 m water depth). 
The central Arctic core, Site B28, was obtained from Lomonosov Ridge (88.87°N, 
140.18°E, 1,990 m water depth). Sediment samples were taken every 1 cm. 
Approximately 1,500 Neogloboquadrina pachyderma sinistral (Ehrenberg) tests 
(5 to 7 mg) were picked from the 150 to 300 um size fraction under a binocular 
microscope. Foraminifera-bound nitrogen isotope analyses were performed at 
Princeton University following the procedures described in ref. (8), which are 
modified from refs. (59) and (60). Briefly, N. pachyderma samples were gently 
crushed; subjected to clay removal, reductive cleaning, and oxidative cleaning; 
dissolved in 4 N HCI to release bound organic nitrogen; oxidized to nitrate with 
basic potassium persulfate; and converted to N,O gas via the denitrifier method 
(61). N isotopes were measured on N,O with a custom-built, automated, helium 
continuous flow-based extraction and purification system coupled to a MAT253 
isotope ratio mass spectrometer (62-64). The analytical precision based on 
long-term replication of internal carbonate-bound organic N standards is 
<+0.30%bo (1sd). The analytical precision of full procedural replicate 5 Nup. 
analyses averaged +0.22%o (1sd, n = 42; error bars plotted on Figs. 2D and 
3A). Further information on sediment core age models, constraints on the iso- 
topic composition of nitrate in the polar Northern Hemisphere (Fig. 1A), and 
the interpretation of Bering Strait sea-level history from 5 Ny, is provided 
in the SI Appendix, Texts S1-S3. 

The growth and decay of ice sheets produces a complex sea-level change 
pattern as the solid Earth responds through crustal deformation, perturbations 
to the Earth's gravitational field, and changes to the Earth's rotation axis. Our 
gravitationally self-consistent sea-level calculations are based on the theory and 
pseudo-spectral algorithm of (65) with a spherical harmonic truncation at degree 
and order 256. These calculations include the impact of load-induced Earth rota- 
tion changes on sea level (66, 67), evolving shorelines, and the migration of 
grounded, marine-based ice (65, 68-70). Our numerical predictions require 
models for Earth’s viscoelastic structure and the history of global ice cover. We 
use an Earth model characterized by a lithospheric thickness of 48 km, and an 
upper and lower mantle viscosity of 0.5 x 10°" Pas and 5 x 10”! Pas, respectively 
(as in ref. (11)). An evaluation of these results using an alternative Earth model 
(VM2) is provided in the S/ Appendix, Text S4. 

Our relative sea-level predictions for the Bering Strait additionally require 
global ice sheet histories to be input into the Earth model described above. 
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Here, we employ three ice sheet histories: ICE-5G, ICE-PC, and W02-derived 
(SI Appendix, Table S1). The ICE-5G ice history is from ref. (19). The ICE-PC ice 
history is from ref. (23) and is constrained such that the LGM and deglaciation 
(26 ka to 0 ka) global mean sea-level history is identical to the ICE-6G ice history 
(71), while the pre-LGM ice history is scaled to fit a global mean sea-level (GMSL) 
history based on fitting GMSL highstand constraints during MIS 3, MIS 5a, and 
MIS 5c (23). The W02-derived ice history is a global ice sheet history for which 
the deglacial ice sheet history is identical to the ICE-5G deglacial history (26 ka 
to 0 ka, ref. 15), while the glaciation phase (before 26 ka) adopts ice geometries 
based on Waelbroeck et al.'s GMSL history (18). Before 26 ka, the ice history at 
each time point is set by the GMSL history in ref. (18), with ice geometry assumed 
to be identical to the post-LGM ICE-5G ice history in the period with the same 
GMSL value. We note that the resultant GMSL history of the W02-derived ice 
history (Fig. 3B) is not identical to the originally published GMSL history of ref. 
(18). The Earth model is run using these three ice histories at 1 to 2 kyr timesteps. 

Because the melt and growth history of the proximal Cordilleran Ice Sheet 
(CIS), and especially its rapid collapse events, can have a first-order effect on 
Bering Strait relative sea level (11), it is also essential to consider possible CIS 
geometries in modeling Bering Strait relative sea level. To address the sensitiv- 
ity of Bering Strait relative sea level to the CIS, we test a suite of CIS ice sheet 
geometries, which incorporate what we term a small, intermediate, and large 
CIS. We produce a set of three ice histories corresponding to the ICE-PC GMSL 
history, which itself is based on the ICE-6G deglaciation history (71) (Fig. 3B and 
Sl Appendix, Table S1). For the intermediate CIS (ICE-PC, Fig. 3C), the ice history's 
geographic distribution assumes that CIS ice geometry prior to 26 ka is identical 
to post-26 ka ice geometry for the same GMSL value based on ICE-6G deglaciation 
history (71). In contrast, for the large CIS (ICE-PC2), prior to 26 ka, the CIS main- 
tains the same boundary as it had during the LGM, with its ice thickness scaled to 
fit the ICE-PC GMSL history (SI Appendix, Fig. S7). This large CIS, characterized by 
maximum ice extent, is part of a North American Ice Sheet history with a reduced 
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MIS 3 eastern Laurentide Ice Sheet (23). For the small CIS (ICE-MIST), we adopt 
the CIS history (defined as west of 120° W, from 80 to 26 ka) from PaleoMIST(24). 
Each of these ice histories is characterized by an identical GMSL history (ICE-PC 
in Fig. 3B) that is maintained by changing ice thickness over Scandinavia and 
Antarctica, locations that are far-field to our site of interest at the Bering Strait. 
For the intermediate and large CIS history (ICE-PC and ICE-PC2), the deglaciation 
history (26 to 0 ka) is adopted from the Gl-31 ice history (11, 30), and incorporates 
substantial melting of the North American ice sheet saddle, which connected the 
Cordilleran and Laurentide Ice Sheets, from 13 to 11.5 ka(S/Appendix, Table S1). 


Data, Materials, and Software Availability. Foraminifera-bound N isotopes, 
sediment core age models, and sea level simulations data are included in the 
SI Appendix, Dataset S1. 


ACKNOWLEDGMENTS. We thank Sergey Oleynik and Laura Gemery for 
laboratory assistance and Peter Huybers for discussions. J.R.F., D.M.S., and 
O.M.U. were supported by NSF OCE-2054780; J.R.F. was also supported by 
the Max Planck Society (with A.M.-G. and G.H.H.) and by the Tuttle Fund of 
the Department of Geosciences, Princeton University (with D.M.S.). T.P. was 
supported by NSF OCE-2054757. T.M.C. was funded by the United States 
Geological Survey Climate Research and Development Program. Any use of 
trade, firm, or product names is for descriptive purposes only and does not 
imply endorsement by the United States Government. 


Author affiliations: “Department of Geosciences, Princeton University, Princeton, NJ 
08544; Max Planck Institute for Chemistry, Mainz 55128, Germany; ‘Earth and Planetary 
Sciences, University of California-Santa Cruz, Santa Cruz, CA 95064; “Department of 
Earth and Planetary Sciences, Harvard University, Cambridge, MA 02138; “Department 
of Marine Sciences, University of Connecticut, Groton, CT 06340; ‘Florence Bascom 
Geoscience Center, United States Geological Survey, Reston, VA 20192; Department 
of Geosciences, Environment and Society, Université Libre de Bruxelles, Brussels 1050, 
Belgium; and "Department of Earth Sciences, ETH Zurich, Zurich 8092, Switzerland 


23. T. Pico, J. R. Creveling, J. X. Mitrovica, Sea-level records from the US mid-Atlantic constrain 
Laurentide ice sheet extend during marine isotope stage 3. Nat. Commun. 8, 15612 (2017). 
24. E.J. Gowan et al., Anew global ice sheet reconstruction for the past 80,000 years. Nat. Commun. 12, 
1199 (2021), 10.1038/s41467-021-21469-w. 
25. P. Fretwell et al., Bedmap2: Improved ice bed, surface and thickness datasets for Antarctica. The 
Cryosphere 7, 375-393 (2013), 10.5194 /tc-7-375-2013. 
26. C.R. Stokes et al., On the reconstruction of palaeo-ice sheets: Recent advances and future 
challenges. Quat. Sci. Rev. 125, 15-49 (2015). 
27. F.D. Hibbert et al., Coral indicators of past sea-level change: A global repository of U-series dated 
benchmarks. Quat. Sci. Rev. 145, 1-56 (2016). 
28. de Gelder, G., Husson, L., Pastier, A-M., Fernandez-Blanco, D., Pico, T., Chauveau, D., Authemayou, 
C., & Pedoja, K. (2022). High interstadial sea levels over the past 420ka from the Huon Peninsula, 
Papua New Guinea. Communications Earth & Environment, 3. https://doi.org/10.1038/s43247-022- 
00583-7 
29. North Greenland Ice Core Project Members, High-resolution record of Northern hemisphere climate 
extending into the last interglacial period. Nature 431, 147-151 (2004). 
30. T. Pico, A. Robel, E. Powell, A. C. Mix, J. X. Mitrovica, Leveraging the Rapid Retreat of the Amundsen 
Gulf Ice Stream 13,000 Years Ago to Reveal Insight Into North American Deglaciation. Geophys. Res. 
Let. 46, 12101-12107 (2019), 10.1029/2019GL084789. 
31. L.E. Lisiecki, J. V. Stern, Regional and global benthic d'°0 stacks for the last glacial cycle. 
Paleoceanography 31, 1368-1394 (2016). 
32. Lehmann, M. F., Sigman, D. M., McCorkle, D. C., Brunelle, B. G., Hoffmann, S., Kienast, M., 
Cane, G., & Clement, J. (2005). Origin of the deep Bering Sea nitrate deficit: Constraints from the 
nitrogen and oxygen isotopic composition of water column nitrate and benthic nitrate fluxes. Global 
Biogeochemical Cycles, 19(4), https://doi.org/10.1029/2005gb002508 
33. J. Granger et al., Coupled nitrification-denitrification in sediment of the eastern Bering Sea shelf 
eads to '°N enrichment of fixed N in shelf waters. J. Geophys. Res. 116, C11006 (2011). 
34. J. Granger, D. M. Sigman, J. G. Gagnon, J.-E. Tremblay, A. Mucci, On the properties of the Arctic 
halocline and deep water masses of the Canada Basin from nitrate isotope ratios. J. Gephys. Res. 
Oceans 123, 5443-5458 (2018). 
35. F. Fripiat et al., Influence of the bordering shelves on nutrient distribution in the Arctic halocline 
inferred from water column nitrate isotopes. Limnol. Oceanogr. 63, 2154-2170 (2018). 
36. Z. W. Brown, K. L. Casciotti, R. S. Pickart, J. H. Swift, K. R. Arrigo, Aspects of the marine nitrogen cycle 
of the Chukchi Sea shelf and Canada Basin. Deep Sea Res. Pt. 118, 73-87 (2015). 
. Ren, D. M. Sigman, R. C. Thunell, M. G. Prokopenko, Nitrogen isotopic composition of planktonic 
oraminifera from the modern ocean and recent sediments. Limnology and Oceanography, 57, 
1011-1024 (2012). 10.4319/lo.2012.57.4.1011. 
38. M.A. Altabet, R. Francois, Sedimentary nitrogen isotopic ratio as a recorder for surface ocean nitrate 
utilization. Glob. Biogeochem. Cy. 8, 103-116 (1994). 
39. A. Randelhoff et al., Pan-Arctic Ocean primary production constrained by turbulent nitrate fluxes. 
Front. Mar. Sci. 7, 150 (2020), 10.3389/fmars.2020.00150. 
J.W. McClelland et al., Particulate organic carbon and nitrogen export from major Arctic rivers. Glob. 
Biogeochem. Cy. 30, 629-643 (2016), 10.1002/2015GB005351. 
hibodeau, D. Bauch, M. Voss, Nitrogen dynamic in Eurasian coastal Arctic ecosystem: Insight 
from nitrogen isotope. Glob. Biogeochem. Cy. 31, 836-849 (2017). 


w 
N 
=. 


41. 


œ% 


pnas.org 


Downloaded from https://www.pnas.org by 178.204.250.13 on January 11, 2023 from IP address 178.204.250.13. 


PNAS 2023 Vol.120 No.1 


M. Straub et al., Changes in North Atlantic nitrogen fixation contro 
501, 200-203 (2013), 10.1038/nature12397. 


lled by ocean circulation. Nature 


B. G. Brunelle et al., Evidence from diatom-bound nitrogen isotopes for subarctic Pacific stratification 


during the last ice age and a link to North Pacific denitrification ch 
PA1215 (2007). 
J. Seguinot, |. Rogozhina, A. P. Stroeven, M. Margold, J. Kleman, N 
Cordilleran ice sheet through the last glacial cycle. The Cryosphere 
K. Lambeck, H. Rouby, A. Purcell, Y. Sun, M. Sambridge, Sea level a 
Last Glacial maximum to the holocene. Proc. Natl. Acad. Sci. U.S.A. 
10.1073/pnas.1411762111. 


anges. Paleoceanography 22, 


umerical simulations of the 
10, 639-664 (2016). 
nd global ice volumes from the 
111, 15296-15303 (2014 


D 


A. S. Dalton et al., Was the Laurentide ice sheet significantly reduced during marine isotope stage 3? 


Geology 47, 111-114 (2019), 10.1130/G45335.1. 


C. L. Batchelor et al., The configuration of Northern hemisphere ice sheets through the quaternary. 


Nat. Commun. 10, 3713 (2019), 10.1038/s41467-019-11601-2. 
P. Liautaud, P. Huybers, Increased sea level sensitivity to CO, forcin 


jcli-d-21-0192.1. 


. J.D. Shakun, D. W. Lea, L. E. Lisiecki, M. E. Raymo, An 800-kyr recor 
and implications for ice volume-temperature coupling. Earth Planet. Sci. Lett. 426, 58-68 (2015). 
M. Suwa, M. L. Bender, Chronology of the Vostok ice core constrained by 0,/N; ratios of occluded air, 
and its implication for the Vostok climate records. Quat. Sci. Rev. 27, 1093-1106 (2008). 


P. C. Tzedakis, M. Crucifix, T. Mitsui, E. W. Wolff, A simple rule to de 


C. W. Snyder, Evolution of global temperature over the past two m 
2016), 10.1038/nature19798. 
F. Parrenin et al., Synchronous change of atmospheric CO, and An 


ransition from ice-albedo and ice-volume nonlinearities. J. Clim. 34, 9693-9708 (2021), 


g across the middle pleistocene 
0.1175/ 


d of global surface ocean 5'0 


ermine which insolation cycles 


ead to interglacials. Nature 542, 427-432 (2017), 10.1038/nature21364. 


illion years. Nature 538, 226-228 


arctic temperature during the last 


deglacial warming. Science 339, 1060-1063 (2013), 10.1126/science.1226368. 


E. Monnin et al., Evidence for substantial accumulation rate variab 


ility in Antarctica during the 


Holocene, through synchronization of CO, in the Taylor Dome, Dome C and DMLice cores. Earth 


Planet. Sci. Lett. 224, 45-54 (2004). 


. T.K. Bauska, S.A. Marcott, E. J. Brook, Abrupt changes in the globa 


glacial period. Nat. Geosci. 14, 91-96 (2021). 


carbon cycle during the last 


P. Köhler, C. Nehrbass-Ahles, J. Schmitt, T. F. Stocker, H. Fischer, A 156 kyr smoothed history of the 


atmospheric greenhouse gases CO, CH,, and N,0 and their radia 
363-387 (2017), 10.5194/essd-9-363-2017. 


e22067421 


ive forcing. Earth Syst. Sci. Data 9, 


19 


57. 


58. 


59. 


60. 


61. 


62. 


63. 


64. 


65. 


66. 


67. 


68. 


69. 


70. 


IA: 


V. V. Pitulko et al., Early human presence in the Arctic: Evidence from 45,000-year-old mammoth 
remains. Science 351, 260-263 (2016), 10.1126/science.aad0554. 

M. Sikora et al., The population history of northeastern Siberia since the Pleistocene. Nature 570, 
182-188 (2019), 10.1038/s41586-019-1279-z. 

A. N. Knapp, D. M. Sigman, F. Lipschultz, N isotopic composition of dissolved organic nitrogen and 
nitrate at the Bermuda Atlantic time-series study site. Global Biogeochem. Cycles 19, GB1018 
(2005). 
H. Ren et al., Foraminiferal isotope evidence of reduced nitrogen fixation in the ice age Atlantic 
Ocean. Science 323, 244-248 (2009). 

D. M. Sigman et al., A bacterial method for the nitrogen isotopic analysis of nitrate in seawater and 
freshwater. Anal. Chem. 73, 4145-4153 (2001). 

K. L. Casciotti, D. M. Sigman, M. G. Hastings, J. K. Böhlke, A. Hilkert, Measurement of the oxygen 
isotopic composition of nitrate and seawater and freshwater using the denitrifier method. Anal. 
Chem. 74, 4905-4912 (2002). 


M. R. Mcllvin, K. L. Casciotti, Technical updates to the bacterial method for nitrate isotopic analyses. 


Anal. Chem. 83, 1850-1856 (2011). 
M. A. Weigand, J. Foriel, B. Barnett, S. Oleynik, D. M. Sigman, Updates to instrumentation and 


protocols for isotopic analysis of nitrate by the denitrifier method. Rapid Commun. Mass Spectrom. 


30, 1365-1383 (2016). 

R.A. Kendall, J. X. Mitrovica, G. A. Milne, On post-glacial sea level - II. Numerical formulation and 
comparative results on spherically symmetric models. Geophys. J. Int. 161, 679-706 (2005), 
10.1111/j.1365-246X.2005.02553.x. 

G. Milne, J. X. Mitrovica, Postglacial sea-level change on a rotating Earth: First results from a 
gravitationally self-consistent sea-level equation. Geophys. J. Int. 126, F13-F20 (1996). 

J.X. Mitrovica, J. Wahr, I. Matsuyama, A. Paulson, The rotational stability of an ice-age earth. 
Geophys. J. Int. 161, 491-506 (2005). 

P. Johnston, The effect of spatially non-uniform water loads on prediction of sea-level change. 
Geophys. J. Int. 114, 615-634 (1993). 

K. Lambeck, A. Purcell, P. Johnston, M. Nakada, Y. Yokoyama, Water-load definition in the glacio- 
hydro-isostatic sea-level equation. Quat. Sci. Rev. 22, 309-318 (2003). 

G.A. Milne, J. X. Mitrovica, J. L. Davis, Near-field hydro-isostasy: The implementation of a revised 
sea-level equation. Geophys. J. Int. 139, 462-482 (1999). 

W. R. Peltier, D. F. Argus, R. Drummond, Space geodesy constrains ice age terminal deglaciation: 
The global ICE-6G_C (VM5a) model. J. Gephys. Res. Solid Earth 120, 450-487 (2015), 
10.1002/2014JB011176. 


https://doi.org/10.1073/pnas.2206742119 


7 of 7 


